file_name = input('请输入图像文件名（包括扩展名）：', 's');


option=input('请选择操作\n1.添加椒盐噪声(噪声概率为0.02)\n2.添加高斯噪声(均值为0.4，标准差为0.01)\n3.滤波算法对比\n4.双边滤波\n');


input_img = imread(file_name);

if option==1
    salt_pepper_noise_img = add_salt_and_pepper_noise(input_img, 0.02); % 噪声概率为0.02
    subplot(1, 2, 1);
    imshow(input_img);
    title('原始图像');
    subplot(1, 2, 2);
    imshow(salt_pepper_noise_img);
    title('添加椒盐噪声后的图像');

elseif option==2

    gaussian_noise_img = add_gaussian_noise(input_img); % 均值为0.4，标准差为0.01
    subplot(1, 2, 1);
    imshow(input_img);
    title('原始图像');
    subplot(1, 2, 2);
    imshow(gaussian_noise_img);
    title('添加高斯噪声后的图像');

elseif option==3

    mean_filtered_image = mean_filter(input_img, 1);
    disp("均值滤波:")
    [PSNR, MSE]=ps(mean_filtered_image,input_img);
    disp(['PSNR: ', num2str(PSNR), ' dB']);
    disp(['MSE: ', num2str(MSE)]);
    % 高斯滤波
    disp("高斯滤波:")
    sigma = 1.5; % 高斯滤波器的标准差
    gaussian_filtered_image = gaussian_filter(input_img, sigma);
    [PSNR, MSE]=ps(gaussian_filtered_image,input_img);
    disp(['PSNR: ', num2str(PSNR), ' dB']);
    disp(['MSE: ', num2str(MSE)]);
    % 中值滤波
    disp("中值滤波:")
    [filtered_img, psnr_value, ssim_value] = median_filter_color(input_img, 3);
    subplot(2, 2, 1);
    imshow(input_img);
    title('原始图像');
    subplot(2, 2, 2);
    imshow(mean_filtered_image);
    title('均值滤波');
    subplot(2, 2, 3);
    imshow(gaussian_filtered_image);
    title('高斯滤波');
    subplot(2, 2, 4);
    imshow(filtered_img);
    title('中值滤波');
elseif option ==4


    X = double(rgb2gray(input_img))/255;
    output_image = Bilater_Gray(X,3,3,0.1);
    subplot(1, 2, 1);
    imshow(input_img);
    title('原始图像');
    subplot(1, 2, 2);
    imshow(output_image);
     title('双边滤波');
end


function output_img = add_salt_and_pepper_noise(input_img, prob)

output_img = input_img;
[height, width] = size(input_img);

% 生成椒盐噪声
for i = 1:height
    for j = 1:width
        if rand() < prob
            output_img(i, j) = 0; % 盐噪声
        elseif rand() < prob
            output_img(i, j) = 255; % 椒噪声
        end
    end
end
end

function output_img = add_gaussian_noise(input_img)

output_img=imnoise(input_img,'gaussian',0.4,0.01);

end


function output_image = mean_filter(input_image, filter_size)

h = fspecial('average', filter_size);
output_image = imfilter(double(input_image), h, 'replicate');
end

function outputImage = gaussian_filter(inputImage, sigma)



outputImage = imgaussfilt(inputImage, 2); % 使用标准差为2的高斯滤波器进行滤波

end

function [filtered_img, psnr_value, ssim_value] = median_filter_color(img, filter_size)

    if size(img, 3) ~= 3
        error('输入图像必须是彩色图像');
    end
    
    % 对每个颜色通道分别进行中值滤波
    filtered_img = img;
    for c = 1:3
        filtered_img(:, :, c) = medfilt2(img(:, :, c), [filter_size, filter_size]);
    end

    % 计算PSNR值
    mse_value = mean((double(img) - double(filtered_img)).^2, 'all');
    psnr_value = 10 * log10(255^2 / mse_value);

    % 计算SSIM值
    [ssim_value, ~] = ssim(filtered_img, img);

    % 显示结果
    fprintf('PSNR值: %.2f dB\n', psnr_value);
    fprintf('SSIM值: %.4f\n', ssim_value);

end
%适用于单通道图像的双边滤波程序

function B = Bilater_Gray(A,w,sigma_d,sigma_r)


[X,Y] = meshgrid(-w:w,-w:w);

G = exp(-(X.^2+Y.^2)/(2*sigma_d^2));


%计算值域核H 并与定义域核G 乘积得到双边权重函数F
dim = size(A);
B = zeros(dim);
for i = 1:dim(1)
   for j = 1:dim(2)

         % 确定作用区域
         iMin = max(i-w,1);
         iMax = min(i+w,dim(1));
         jMin = max(j-w,1);
         jMax = min(j+w,dim(2));
         %定义当前核所作用的区域为(iMin:iMax,jMin:jMax)
         I = A(iMin:iMax,jMin:jMax);%提取该区域的源图像值赋给I

         %计算值域核H.
         H = exp(-(I-A(i,j)).^2/(2*sigma_r^2));

         % Calculate bilateral filter response.
         F = H.*G((iMin:iMax)-i+w+1,(jMin:jMax)-j+w+1);

         B(i,j) = sum(F(:).*I(:))/sum(F(:));

   end
end

end
